Generalized method of image and the tunneling spectroscopy in high-T c 

superconductors 



(N 
O 

o 



in 



Shin-Tza Wu 1 and Chung- Yu Mou 1 ' 2 
1. Department of Physics, National Tsing Hua University, Hsinchu 30043, Taiwan 
2. National Center for Theoretical Sciences, P.O.Box 2-131, Hsinchu, Taiwan 

(Dated: February 1, 2008) 

A generalized method of image is developed to investigate the tunneling spectrum from the metal 
into a class of states, with the tight-binding dispersion fully included. The broken reflection sym- 
metry is shown to be the necessary condition for the appearance of the zero-bias conductance peak 
(ZBCP). Applying this method to the d-wave superconductor yields results in agreement with ex- 
periments regarding the splitting of ZBCPs in magnetic field. Furthermore, a ZBCP is predicted 
for tunneling into the (110) direction of the d-density wave state, providing a signature to look for 
in experiments. 
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The current transport through a heterojunction con- 
sisting of a normal metal and another different material 
(X) has been the subject of interest for many years. In 
this setup, the normal metal with known spectral prop- 
erties is used as a probe to analyze the electronic states 
of the material X.a Although such measurements have 
provided useful insights into the bulk spectral proper- 
ties of X, it has been also realized that the presence of 
the interface matters. The zero-bias conductance peak 
(ZBCP) observed in the tunneling spectra when X is a 
d- wave superconductor (ND junction) in (110) directiono 
is a well-known example of interface effects. However, 
the issue of exactly how the tunneling measurements are 
related to the bulk properties has never been answered 
satisfactorily. Conventionally, the ND junction is ana- 
lyzed in the mean-field level, using the Bogoliubov de 
Gennes (BdG) equations in which continuum and quasi- 
classical approximations are often invoked. While these 
approximations are valid for conventional superconduc- 
tors, they are certainly not justified for high Tc cuprates 
where proximity to the Mott insulators entails fully, con- 
sideration of the tight-binding nature. PreviouslytJ this 
was done by numerically solving the discrete BdG equa- 
tion for each interface orientation individually without 
elucidating their relations to the bulk properties. This 
technical inconvenience makes it difficult to include fluc- 
tuations systematically in this approach. 

In this work, we shall adopt a different approach based 
on the non-equilibrium Keldysh-Green's function formal- 
ism which enables one to construct systematically higher 
order cor*fi&fcions from the mean field lattice Green's 
functions D'EIeI In this approach, because X extends over a 
semi-infinite space, one shall need the half-space Green's 
functions. For simple configurations such as the (100) 
orientation of a d-wave superconductor, it turns out that 
this half-space Green's functions only differ from the bulk 
ones by sinusoidal factors. This relation certainly does 
not hold for other orientations as it predicts no ZBCP 
in the (110) direction. We shall develop a generalized 
method of image which enables us to construct the half- 
space Green's function from the bulk ones. We empha- 



size the generality of this method and its ability to account 
for the low energy features in the tunneling spectrum for 
a whole class of states. As a demonstration, in this ar- 
ticle we will focus mostly on the study of ND junctions 
and only briefly mention the applications to other sys- 
tems. The effects of interactions and fluctuations will be 
addressed elsewhere. 

Our results indicate broken reflection symmetry is nec- 
essary for the emergence of ZBCP. For ND junctions our 
method can reproduce earlier results on the ZBCP in the 
continuous wave approximation.oQ In a full tight-binding 
calculation for (110) and (210) directions, we obtain the 
doping dependence of the ZBCPs which shows its sensi- 
tivity to the Fermi surface topology. In particular, the 
splitting of the ZBCP in the current-carrying state is also 
calculated and is shown to be in agreement with exper- 
iments. At the end, we analyze the case when X is the 
d-density wave (DDW) state in (110) direction and the 
semi-infinite graphene sheet with a zig-zag type interface. 
The former state was recently pepposed as a possible nor- 
mal state for high Tc cuprates E Conductance peaks are 
found for both states. 

We start by considering a junction consisting of a 2D 
normal metal on the left [L) hand side (— oo < x < —d, 
d is the lattice constant of the metal side) and a e?-wave 
superconductor (0 < x < oo) on the right (R) hand side 
(see Fig. 1), governed by the Hamiltonian Hl and Hr 
respectively. At the mean-field level, we have 
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where < ij > denotes the nearest-neighbor (NN) bond, 
< ij >' the next NN bond, and Ay possesses the d-wave 
symmetry. The tunneling Hamiltonian connects the in- 
terface points at x = — d and x = 0, and is given by 

Ht = J2 y - Vr\)(cl c R + CrCl), where the summa- 
tion is over lattice points along the interface, chosen to 
be in the y direction. We shall assume that both sides are 
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FIG. 1: Schematic plot of the ND junction in (110) and (210) 
directions. Here two hard walls for the (210) case are at A = 
-2d/y/E and B = -2d/VE. The dots on the x axis are the 
reduced ID lattice points. 

square lattices and characterize the orientation of RHS by 
the Miller indices (hkO). The total grand Hamiltonian is 
then given by K = (H L - fi L N L ) + (H R - fi R N R ) + H T = 
Kq+Ht- Here and ^i R are the chemical potentials and 
their difference /i^ — fi R is fixed to be the voltage drop 
eV across the junction . 

In the non-equilibrium Keldysh-GreenJa. function for- 
malism, Ht is adiabatically turned on.Dtj As a result, 
the bare Green's function is denned only on a half plane. 
Since the nearest-neighbor bonds to the interface sites are 
cut, there is effectively a hard wall located at say, for the 
(110) interface, x = —d/2. This hard-wall boundary con- 
dition prompts the application of the method of image. 
However, because lattice points in the half-plane usually 
do not form a simple Bravais lattice and the d-wave gap 
changes sign under reflection, the implementation of the 
conventional method of image appears problematic. To 
overcome these difficulties, a Fourier transform in the y 
direction is performed first. Consider the case of (110) 
orientation, the Hamiltonian with only NN hopping be- 
comes 

Hr = - 2t R C0S (^Yj 4*( k y) c *+i<?( k y) 

+ai(-k y )a+ ir (ky)} + R.c., (2) 

where —n/d < k y < ir/d and 2Ao is the gap value. 
The whole problem is now one dimensional, and the 
hard wall becomes a point. Note that the suppression 
of the gap near the interface can be taken into account 



by adding self-consistent 6A Q (i) to Eq. (0), and can 
be treated perturbatively later. In the presence of t' R , 
additional terms J2i,a,k y -2t' R cos(k y d) c\ a {k v )c la {k v ) + 

t'it c \ ( y{ k y) c i+2a{k y ) + H.c. appear. Since at the bound- 
ary, both t R and t' R are cut and become dangling bonds, 
one needs to introduce two hard walls at x = —d/2 and 
x = — d. For clear presentation, we shall first set t' R = 0. 
In this case, we are looking for the Green's function, 
k y x, x') (which is a 2 x 2 matrix in Nambu's no- 
tations), that satisfies the boundary condition Go = at 
x = —d/2. We shall suppress the dependence on u> and 
k y . Here since x' is the location of the point source and 
its image point is at — d — x'p the method of image can 
be employed by constructing 

G (x,x') = G (x - x') - G (x + d + x') a(x'), (3) 

where Go is the bulk bare Green's function and a is a ma- 
trix to be determined. The first term is the direct propa- 
gation from x' to x, while the second term will reduce to 
the propagation from the image point to x in special situ- 
ations (see below). In fact, since Go has to vanish at x — 
-d/2, we obtain a(x') = G^\d/2 + x')G (-d/2 - x'). 
Therefore, the second term describes the propagation 
from x' to x via the reflection of the hard wall. The 
matrix a, apart from fitting the boundary condition, car- 
ries the important information about the gap structure 
along the reflected path from x' to x. Note that in cal- 
culating the tunneling current, since Ht only connects 
points along the interface, only the surface Gjjpep's func- 
tion g (u!,k y ) = Gq(x = 0, x' = 0) is neededlM Writing 
Go in the Fourier k x space, we find 

f 27T / d dk x 

go(uj,k y )= / —^G (uj,k y k x )[l-exp(ik x d)ao}, 
J-27r/d 27r 

(4) 

where the factor otQ = Gq 1 (d/2)Go(— d/2) is indepen- 
dent of k x . If the reflection symmetry holds for the state 
X, such as d-wave superconductors in (100) direction 
(in this case, d/2 is replaced by d), onft-,has a = 1 
and Eq. (||) reduces to the familiar formQ g~o(uj,ky) = 

/-27r/d dk x /ii Gq(uj, k y k x ) sin 2 (k x d/2). Therefore, apart 
from modifications due to the sinusoidal factors, the den- 
sity of state almost has the same feature as the bulk one. 
However, for other orientations such as the (110) direc- 
tion, reflection symmetry with respect to the interface is 
broken. As a result, ao is not the identity matrix and as 
we shall see, this will give rise to the ZBCP. 

The advantage of Eq. (^) is that it is purely based 
on the bulk Green's functions. The interface orientation 
is encoded in k x and k y . In other words, Go(oj,k y ,k x ), 
which appears in Eq. (ji|) and «o, is simply the usual 
bulk BCS Green's function but with k being rotated by 
45°. This technical merit is retained for other interface 
orientations with k being rotated by the angle in accor- 
dance with the interface orientation. More importantly, 
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this also offers a scheme for studying fluctuations and in- 
teractions. Essentially one can take into account these 
effects through the bulk Green's function Go- This will 
be explained in more detail in a separate publication. 

When evaluating Gq(x — x'), the dominant contribu- 
tions come from the poles determined by (aj + irj) 2 — E 2 = 
0, where Ek = \J('k + A 2 ,. In the continuum limit, the 
dispersion becomes e& = h 2 (k 2 — k 2 Fx )/2m and A& = 
A 0082(0-6*0), where k 2 Fx = k 2 F ~k 2 , 9 = sin~ 1 (k y /k F ) 1 
and 6*o is the angle between the crystal a axis and x direc- 
tion. At the same time, the integration range of k x is ex- 
tended to ±oo. There are four poles located at ±fc± with 

k± = \J k Fx ± 2my / (w + it]) 2 — \A±\ 2 /h 2 , representing 
particles and holes along different directions. Here A± 
are gaps in directions ±0. By contour integration, one ob- 
tains Gq(x — x') and thuSpQ^x'). After some algebra and 
assuming that k F is large, til indeed Eq.(||) reproduces re- 
sults obtained in Ref . by direct solving the equations of 
motion. In our approach, the continuum approximation 
is not necessary. To investigate any effect that is due to 
the tight binding nature, the full tight binding disper- 
sion has to be retained. In this case, the integration over 
k x can not be extended to ±oo, and thus poles are in 
different structure and a substantial difference from the 
continuum approximation could result. 

We now include the hopping t' R for the (110) direction. 
The main complication is to add a second hard wall at 
x = —d. This is a simple generalization of the single 
hard-wall problem. One simply requires Go vanish on all 
these hard walls simultaneously. Therefore, we write 

Gq(x, x) — Gq(x — x') — Gq{x — x' 1 )ai(x') 

- G (x - x' 2 )a 2 (x'), (5) 

where x[ = —x 1 — d and x' 2 — —x' — 2d are images of 
x' . The boundary conditions at x = —d and x = —2d 
determines a\ and a 2 . The surface Green's function thus 
obtained is the bare one and will get renormalized by Ht, 
giving rise to-feur different components in the differential 
conductanceolill. The strength of Ht, characterized by 
t, determines the relative weight among each component. 
In Fig. 2, we show our results for the spectrum of the to- 
tal differential conductance for various dopings. The pa- 
rameters adopted are determined self-consistently from 
the mean-field slave-boson theory for the t-t'-J model.cl 
It is clear that the ZBjCPs are the most important fea- 
tures at low energies.u Since the ZBCP arises from the 
existence of zero-energy states, it must appear as poles 
at zero energy in the Green's function. For (110) direc- 
tion without t' R , this is entirely determined by zeros of 
the denominator in a: 

P(u J ,k y ) = det[Go(d/2)}. (6) 

In the continuum limit (d — ► 0), (3 can be evaluated ana- 
lytically: (3 « -[cj 2 /A 2 + (A+/|A+| + A_/|A_|) 2 ] with 
A = |A+||A_|/(|A+| + |A_|). Therefore, poles of g 
at lu = depends crucially on whether there is a sign 



0.14 




-0.5 -0.2 0.0 0.2 0.5 

V (in unit of J) 



FIG. 2: The total differential conductance of several dop- 
ings for (110) interface with 77 = 0.01 and th = 1.0. 
The weak link is modeled by the interface hopping t(ui) = 
exp(— y/ (uo — |w|)/r) with uiq = llAo and Y = An. Inset: 
dl/dV curve for (210) interface with 5 = 0.08. 



change of the gap on the Fermi surface. This criterion, 
however, does not hold as one goes to the tight binding 
limit because pairing no longer only occurs on the Fermi 
surface. As a result, a numerical computation of (3 is nec- 
essary. Our results indicate that ZBCPs are sensitive to 
the Fermi surface topology. In fact, for (110) surface, the 
height of the peak depends on the volume of the Fermi 
surface. It reaches maximum when (1r — and decreases 
when [ir ^ 0. For other orientations, the ZBCP could 
even disappear. 

For general orientations, there could be more than one 
hard wall. As a demonstration, we consider the (210) 
interface. In this case, when t' R = 0, two hard walls are 
located at x = —dj\fh and — 2d/v / 5, in analogy to the 
(110) surface with t' R . A typical result for small scale of 
t is shown in the inset of Fig. 2. As t increases, the zero- 
energy states are able to leak out, and thus the ZBCPs 
get broadened. Note that lattice points with dangling 
bonds form a pair-breaking region near the interface, re- 
sulting in the peaks around 2Ao- They are due to quasi- 
particle bound states with non-zero energy.113 

The ZBCPs split in the presence of magnetic fields 
H, essentially due to the Doppier shift caused by the 
supercurrent near the interfaceIL3 In the tight-binding 
model, Ajj is shifted to Ay exp[iq • (r. + r^)], where 
q = eHX/2hc with A being the penetration depth. By 
redefining c- la — ci a exp(iq • r^), the dependence on q can 
be absorbed into e&. Figure 3 shows the field dependence 
of the splitting for the (j-10) interface, in agreement to 
recent experimental dataES It is seen that for large H, 
the splitting deviates from linear dependence on H due 
to the lattice effect in our approach. The inset shows the 
doping dependence of splitting, reflecting its sensitive de- 
pendence on the Fermi surface topology. In fact, in the 
special case when particle-hole symmetry holds (for ex- 
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FIG. 3: The field dependence of the splitting for an under- 
doped case (8 = 0.12). Inset: Doping dependence of the 
splitting for a fixed magnetic field. 
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ample, (1r = and t' R = 0), we find that no matter how 
large H is, the ZBCP does not split at all, in consistent 
with naive expect at ions .til 



Finally, our approach is easily modified to deal with 
other states. When X is the DDW, the two-component 
indices are associated with the two sublattices in this 
state. The formulation presented for the,ND case can 
then be applied with minor modification]^ At (110) ori- 
entation, the DDW state does not possess reflection sym- 
metry; one thus expects ZBCP's in this case. Figure 4 
shows a typical result for tunneling into the DDW state in 
the (110) direction when = 0. invariabbz_the ZBCP 
is present, consistent with a recent report .E_3 For finite 
fiR, unlike the <i-wave superconductor, one simply adds 
/ifl to the quasi-particle energy -Efcjj resulting in shifting 
of the ZBCP to the bias at ^lr. The existence of this 
conductance peak thus provides a signature of the DDW 
state. A similar analysis can be done for the semi-infinite 
graphene sheet. When projected onto ID lattices, it is 
then obvious that reflection symmetry is preserved in the 
case of the armchair interface while not in the zigzag case. 
This results in for the latter a ZBCP in the dl/dV curve 
(see the inset of Fig._4)(Ref. |ll|) - consistent with previ- 
ous numerical work.li3 In conclusion, we have developed 
a versatile formulation for studying the tunneling spec- 
tra of junction systems. In particular, it is capable of 
predicting whether the ZBCP shall arise or not in the 
conductance curve. The results discussed in this article 
represent some typical examples; further applications to 
other systems will be reported elsewhere. 



FIG. 4: A typical dl/dV curve for tunneling into the (110) 
direction of a d-density wave state. Here [ir = and Ao = 
0.1. Inset: A similar plot for tunneling from a wide-band 
metal into the graphene sheet with zig-zag interface. Here 
the hopping amplitude 70 =0.1 
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